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Abstract 

In laboratory experiments that impose shear deformation on partially molten 
aggregates of initially uniform porosity, melt segregates into high-porosity sheets 
(bands in cross-section). The bands emerge at 15-20° to the shear plane. A 
model of viscous anisotropy can explain these low angles whereas previous, simpler 
models have failed to do so. The anisotropic model is complex, however, and 
the reason that it produces low-angle bands has not been understood. Here we 
show that there are two mechanisms: (i) suppression of the well-known tensile 
instability, and (ii) creation of a new, shear-driven instability. We elucidate 
these mechanisms using linearised stability analysis in a coordinate system that 
is aligned with the perturbations. We consider the general case of anisotropy that 
varies dynamically with deviatoric stress, but approach it by first considering 
uniform anisotropy that is imposed a priori and showing the difference between 
static and dynamic cases. We extend the model of viscous anisotropy to include 
a strengthening in the direction of maximum compressive stress. Our results 
support the hypothesis that viscous anisotropy is the cause of low band-angles in 
experiments. 


1 Introduction 

In laboratory experiments, forced shear deformation of nominally uniform, partially 
molten rocks causes melt segregation into high-porosity bands oriented at low angle 
(15-20°) to the shear plane [Holtzman et ah, 2003, Holtzman and Kohlstedt, 2007, 
King et ah, 2010]. Stevenson [1989] predicted the emergence of such bands in a self¬ 
reinforcing feedback arising from the porosity-weakening of the crystal-l-magma aggre¬ 
gate, but the angle predicted by this theory was 45° [Spiegelman, 2003], much higher 
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than observed. The low angle of high-porosity bands is widely thought to provide an 
additional constraint on the rheology of the aggregate, but it has proven challenging to 
understand. Katz et ah [2006] found that non-Newtonian viscosity with a high sensitiv¬ 
ity to stress could reproduce the low angle of bands, but King et ah [2010] subsequently 
showed that the viscosity of experiments that produce low-angle bands is actually close 
to Newtonian. 

Theory by Takei and Holtzman [2009a,b] of anisotropic viscosity under diffusion 
creep of a partially molten aggregate represents a possible solution. This theory is 
motivated by observations of the coherent alignment of melt-pockets between solid 
grains under a deviatoric stress [e.g. Daines and Kohlstedt, 1997, Zimmerman et ah, 
1999] and of the enhancement of diffusion creep by melt at grain boundaries and triple 
junctions [e.g. Cooper et ah, 1989]. The melt is a fast pathways for diffusional transport 
of solid constituents around grains; the alignment of melt with respect to the principal- 
stress directions hypothetically results in anisotropic viscosity of the aggregate [Takei 
and Holtzman, 2009a]. 

Analysis of the theory of anisotropic viscosity by Takei and Holtzman [2009b], Takei 
and Katz [2013], Katz and Takei [2013], and Allwright and Katz [2014] shows that it 
introduces qualitatively different behaviour from previous models with isotropic (and 
even power-law) viscosity. Shear and normal components of stress and strain-rate are 
coupled under viscous anisotropy; as a result of this coupling, a gradient in shear stress 
becomes a driving force for melt segregation that is not present in the isotropic system. 
Under Poiseuille flow, melt segregates toward higher-stress regions; under torsional 
flow, compressive hoop stresses drive the solid outward and the magma inward. The 
mechanics of this “base-state” melt segregation are explained in detail by Takei and 
Katz [2013]. An experimental test of radial melt segregation in torsional flow by Qi 
et ah [2015] shows striking consistency with predictions. 

Furthermore, theoretical work has demonstrated that there is a connection between 
the strength of anisotropy and the angle of high-porosity bands that emerge by unstable 
growth. This was shown with linearised stability analysis [Takei and Holtzman, 2009b, 
Takei and Katz, 2013] and numerical simulations [Butler, 2012, Katz and Takei, 2013] 
where the strength and orientation of anisotropy are assumed to be known and are 
imposed a priori. In those static-anisotropy calculations, high-porosity bands emerge at 
low angles to the shear plane only when viscous anisotropy is at or near saturation. This 
is a rather restrictive condition that may be incompatible with the robust appearance 
and consistently low angle of bands in experiments [Holtzman and Kohlstedt, 2007]. 
However, in numerical simulations that allow anisotropy strength and direction to vary 
dynamically in space and time [Katz and Takei, 2013], band angles are signihcantly 
lowered and appear to be less sensitive to the mean strength of anisotropy. These 
hndings raise several basic, unanswered questions: Why do the mechanics of viscous 
anisotropy give rise to low-angle bands? Why is dynamic anisotropy more effective 
in this regard than static anisotropy? What are the general conditions under which 
low-angle, high-porosity bands should form? 

The present manuscript addresses these questions through a combination of lin- 
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earised stability analysis and physical reasoning. The crucial, enabling advance is to 
perform the analysis in a coordinate system that is rotated to align with the porosity 
bands (rather than with the plane of shear). This drastically simplihes the expressions 
for growth rate under static anisotropy [Takei and Katz, 2013], making them readily 
interpretable in physical terms. Moreover, it allows us to extend the analysis to dy¬ 
namic anisotropy in a form that exposes the physical differences from static anisotropy. 
Finally, the same coordinate rotation clarihes the physical reason for low angles under 
isotropic, non-Newtonian viscosity. 

The manuscript is organised as follows. In the next section, we briefly discuss the 
nondimensionalised governing equations and present an anisotropic, viscous constitu¬ 
tive model for the two-phase, partially molten aggregate. The full, non-linear system 
is solved numerically in §3 for static and dynamic cases, to elucidate the questions 
listed above. The coordinate rotation is introduced and the linearised stability analysis 
is developed in §4. In particular, §4.3 develops an expression for the growth-rate of 
porosity perturbations under the fully dynamic model of §2. This expression is chal¬ 
lenging to understand and so we subsequently consider it under reducing assumptions 
of static anisotropy (§5.1), which includes the simplest case of Newtonian, isotropic 
model. We build on this to explain the full complexity in §5.2 and §5.3. We conclude 
with a summary and discussion of the results in terms of the motivating questions. 


2 Governing and constitutive equations 


In the theory of magma/mantle interaction, the macroscopic behaviour of a two-phase 
aggregate is treated within the framework of continuum mechanics [e.g. Drew, 1983, 
McKenzie, 1984]. This theory is concerned with the evolution of macroscopic helds 
including the volume fraction of melt or porosity 0, the velocity of the solid phase V, 
the liquid pressure P (compression positive), and the bulk or phase-averaged stress 
tensor aij = (1 — ~ 4’P^ij where af^ is the stress tensor of the solid phase (tension 

positive). Further details of the two-phase-flow theory were previously presented [e.g. 
Takei and Katz, 2013, Rudge et ah, 2011] and are not repeated here. 

We proceed directly to the nondimensional governing equations. 
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and refer the reader to Takei and Katz [2013] and references therein for details of the 
derivation and rescaling. In the system (1) we have introduced the differential stress 
tensor Tij = aij + PSij. Also, we have excluded body forces and assumed that the 
permeability of the solid matrix is a function of the porosity only, proportional to 
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{.4’14’oY-i where 00 is a reference porosity and £ is a constant. R is the nondimensional 
compaction length and is a rheological parameter explained below. To close the 
system, a constitntive relationship that relates the differential stress Tij and strain rate 
is reqnired. 

Takei and Holtzman [2009b] and Takei and Katz [2013] proposed a model of anisotropic 
viscosity cansed by stress-indnced microstrnctnral anisotropy. In partially molten rocks, 
the melt phase is contained within a permeable network of tnbules between grains. The 
solid matrix is formed by a contignons skeleton of solid grains. The area of grain-to- 
grain contact is known as the contignity. Contignity is the microstrnctnral variable 
that determines the macroscopic (i.e., continnnm) mechanical properties of the matrix 
[Takei, 1998, Takei and Holtzman, 2009a]. Althongh the eqnilibrinm microstrnctnre 
developed nnder hydrostatic stress has isotropic contignity, deviations from the eqnilib- 
rinm microstrnctnre have been observed in experimentally deformed, partially molten 
samples [e.g. Daines and Kohlstedt, 1997, Takei, 2010]. Based on these observations, 
we infer that nnder a differential stress, the grain-to-grain contacts with normals that 
are parallel to the maximnm tensile stress (ts) are rednced in area; similarly, the areas 
of those with normals parallel to the maximnm compressive stress (ri) are increased. 
Using a microstrnctnre-based model of aggregate viscosity [Takei and Holtzman, 2009a] 
and a coordinate transformation [Takei and Katz, 2013], the constitntive law and the 
viscosity tensor are 

Rj ^ijkl^kli 





ci 

kl^ XX 

YY 


XX 

( ^s + |-^cos20 


sin 20 \ 


“/(3 + COS40) 

“/(I cos 40) 

sin 40 

YY 


rg + | + ^cos20 

sin 20 



“/(3 + COS40) 

+ "/ sin 40 

XY 

V 


1- ^^^(1-cos 40) / 


For simplicity, we consider a two dimensional problem, in which the Ti-Ts plane is 
parallel to the X-Y plane. Therefore only the two-dimensional version of Cijki is 
written in (2b). Only 6 of the 16 components are shown dne to the symmetry of Cijki 
nnder the exchange of i and j, k and /, and ij and kl. 

The factor in front of the matrix represents the normalised shear viscosity ?7(0)/?7(0o); 
it decreases exponentially with increasing melt fraction 0 , and so A is called the porosity¬ 
weakening factor. We take A = 27 based on the experimental resnlts [e.g. Mei et al., 
2002]. The parameter represents the bnlk-to-shear viscosity ratio, = ^/rj, which 
is assnmed to be constant (=5/3) based on theoretical results by Takei and Holtzman 
[2009a] [although see Simpson et ah, 2010a,b]. Parameters a, 0 , and 0 represent the 
magnitude and direction of microstrnctnral anisotropy: a and 0 quantify the amplitude 
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of contiguity reduction and increase, respectively, in the directions of principal stress 
Ta and ri; 0 represents the angle that the most-tensile stress (ra) direction makes with 
the X-axis of the coordinate system. Using the local differential stress Tij, 0 is given 
by 

tan 20 =-——, (3) 

Txx — tyy 

and a and {3 are modelled as 

2(Xt Toffset) 

"^sat 

(3 = r^a, (4b) 

where At = ra — ti = \/ {txx — tyyY + represents the amplitude of deviatoric 
stress. The detailed forms of the functions in (4) are poorly constrained, owing to a 
lack of experimental data. The form of a was chosen based on the constraints that a is 
less than or equal to 2 [Takei and Katz, 2013] and increases with increasing differential 
stress [Daines and Kohlstedt, 1997, Takei, 2010]. The parameter is assumed to be a 
constant that is probably between 0 and 1. In the present study, a is parameterised by 
"roffset and Tsat, which control the stress-offset and slope of increase. In Takei and Katz 
[2013], (3 = 0 and a was parameterised by Tgat alone. Parameters (3 and Toffset are newly 
introduced here. 

For simplicity in previous linearised analyses [Takei and Katz, 2013, Allwright and 
Katz, 2014], parameters a, (3, and 0 were hxed to their initial values. We call this 
simplifying assumption the static anisotropy model. In contrast, the complete model 
with stress-dependent direction and magnitude is called the dynamic anisotropy model. 
Katz and Takei [2013] discovered a remarkable difference between static and dynamic 
anisotropy; this difference motivates the present study and is demonstrated in the next 
section. 



a = 1 + tanh 


3 Numerical solutions 

Numerical solutions of equations (1) and (2) highlight the difference between results 
obtained for static and dynamic anisotropy. The solutions are computed with a hnite- 
volume method on a fully staggered grid that is periodic in the X-direction; in this 
section, the X axis is taken parallel to the initial flow direction. No-slip, impermeable 
boundary conditions enforce a constant displacement rate of plus or minus on the 
top and bottom boundaries, respectively. A semi-implicit, Crank-Nicolson scheme is 
used to discretise time and the hyperbolic equation for porosity evolution is solved 
separately from the elliptic system in a Picard loop with two iterations at each time- 
step. The solutions are obtained in the context of the Portable, Extensible Toolkit for 
Scientihc Computation [PETSc, Balay et ah, 2001, 2004, Katz et ah, 2007]. Full details 
and references are provided by Katz and Takei [2013]. 


5 






Figure 1: Comparison of numerical solutions to (1) and (2) with static and dynamic 
anisotropy. For both calculations, = 1, i? = 1, = 5/3, cpo = 0.05, e|0i(X)| < 0.005 

and the domain is discretised into 600x300 square cells, (a) Porosity field at a strain of 
0.75 for a simulation with a = 1.8 and 0 = 45° throughout the domain, (b) Porosity field 
at a strain of 1.25 for a simulation with anisotropy calculated according to (3) and (4) with 
"Tsat = Ij TofTset = 1-5. (c) Spectral power binned by wavefront-angle 6 to the shear plane 
[after Katz et ah, 2006] for the porosity fields shown in panels (a) and (b). Each spectrum is 
normalised by its maximum power, (d) and (e) Two-dimensional histograms derived from 
the simulation with dynamic anisotropy at a strain of 1.25 [after Katz and Takei, 2013, 
Fig. 12]. Red dashed lines have a slope given by the ratio of perturbation quantities 
and Qi/4>i from the stability analysis in §4. 
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Figure 1 compares solutions with fixed and dynamic anisotropy. In panel (a), 
anisotropy parameters are prescribed as a = = 1.8, 0 = 45°; in panel (b), these 

parameters are calculated cell-wise using equations (3) and (4), with = 1. Both 
calculations have R = 1 (compaction length equal to domain height) and are initialised 
with the same porosity field, 0(X, t = 0) = (j)o+e(j)i{X), where 0o = 0.05 and e = 0.005. 
4>i{X) is a smooth, random field with unit amplitude, generated by filtering grid-scale 
white noise to remove variation at wavelengths below 15 grid-cells. Because the growth- 
rate of porosity perturbations differs for fixed and dynamic anisotropy, the simulations 
are shown at different values of the average simple-shear strain 7 . 

The different orientation of high-porosity features is evident in panels (a) and (b): 
dynamic anisotropy is associated with lower angles. This is quantified by the power 
spectrum in panel (c), where the power from a 2D fast-Fourier transform of the porosity 
field is binned according to the angle between the wavefront and the shear plane [Katz 
et ah, 2006]. Dynamic anisotropy produces a peak at ~ 10° whereas static anisotropy 
produces a peak at ~ 23°. There is also a high-angle (~ 80°) peak for static anisotropy 
(corresponding to features visible in panel (a)) that does not survive at large strain. 
Panels (d) and (e) show the covariation of a and 0 with 0 in panel (b); black dotted 
lines indicate mean values. These means are closely matched with the parameter values 
used in the fixed-anisotropy simulation. It is therefore clear that the difference in 
the dominant band angle (panel (c)) arises from the coupling between stress and the 
variations in a and 0. What is unclear, however, is the physical explanation for this 
difference and, indeed, why viscous anisotropy gives rise to bands at angles less than 
45° to the shear plane at all. We clarify these points below. 

4 Linearised analysis with a perturbation-oriented 
coordinate system 

Let e <C 00 be the initial amplitude of a porosity perturbation. We express the problem 
variables as a series expansion about the base-state in which the porosity is uniform 
and equal to 0o. We truncate the series after the first-order terms. 
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The first term of (5) with index 0 represents a simple-shear flow and its associated 
anisotropy, which is the base-state solution of order one (e°), corresponding to the 
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Figure 2: Schematic diagrams of the coordinate axes and porosity perturbation, (a) The 
coordinate system (X, Y) for the linearised analysis is taken such that the X axis is parallel 
to the initial perturbation wavefronts. The shear plane of the base-state, simple-shear flow 
is then rotated by an angle 6. (b) The base-state normal stress Ty^y is oriented parallel to 
the initial perturbation wave vector, (c) The base-state shear stress rj^y is parallel to the 
wavefronts. 


uniform porosity ^q. The second term of (5) with index 1 represents the perturbation 
of order caused by e0i. By substituting (5) into equations (1), using V- = 0, 
and balancing terms at the order of e^, we derive the governing equations for the 
perturbations as 
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(6b) 

(6c) 


where D0i/Df = d(j)i/dt -|- • V0i. 

Following previous studies, the porosity perturbations take the form of a plane 
wave oriented at a given angle to the base-state shear plane. Past workers chose to align 
the coordinate system with the base-state shear plane, such that the base-state strain 
rate tensor has a simple form [e.g. Spiegelman, 2003]. Although the coordinate system 
was so aligned in the numerical models above, in this section the coordinates are rotated 
such that the Y axis is parallel to the wave-vector of the initial perturbation, as shown 
in Figure 2. With this choice, 6 again represents the angle between the perturbation 
wavefronts and the base-state shear plane. However, in the rotated coordinate system, 
0 depends on both the direction of and the orientation of the bands. The base-state 
direction of maximum tensile stress Ts makes an angle vr/d to the shear plane (eqn. (3)) 
and so for a coordinate rotation by 6, we have ©o = vr/4 -|- 6 (Figure 2a). 

In the following part of this section, we give an outline of the linearised approach, 
which shows that the new coordinate system reduces the complexity of the analysis and 
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exposes the physical mechanisms of perturbation growth. This enables us to clarify the 
mechanics leading to low-angle bands in complicated problems such as under dynamic 
anisotropy. We hrst consider the base-state, simple shear flow at the order of e° (§4.1) 
and then the linearised governing equations at the order of (§4.2). Finally, in §4.3, 
we obtain the growth rate of porosity perturbations 01 for the most general case of 
dynamic anisotropy. The result obtained is used in §5 to clarify the mechanisms of 
low-angle-band formation. 


4.1 Base-state simple shear flow 

Using the angle 6 between the initial perturbation wavefronts (aligned with the X 
direction) and the base-state shear plane (Fig. 2a), the components of the base-state 
strain rate tensor in the rotated coordinate system are 

.( 0 ) _ 1 / — sin26' cos26' \ . . 

2 0 cos26' sin 20 y 


As shown in Figs. 2b and 2c, r^y and rj^y represent, respectively, the base-state 
tensile and shear stresses normal and parallel to the perturbation wavefronts, which 
play important roles in understanding the growth of these perturbations. Noting that 
Cyy = —these components are given by 
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( 8 ) 


with 00 = f -|- 0. In Figure 3, Tyy and r^y are plotted as a function of angle 0. An 
understanding of their systematics is needed to interpret the results of the stability 
analysis. 

Considering hrst the solid curves representing normal stress, we see that for = 
00 = 0 (isotropic, top-left panel), Tyy follows the expected pattern of sin 20. It is tensile 
for 0 < 90° and compressive for 0 > 90° [Spiegelman, 2003]. However, as ao increases 
(for 00 = 0, top row of panels), Tyy becomes negative (compressive) at all angles. The 
mechanism for this change is two-fold. First, decreases the viscous resistance to 
extension in the direction and reduces the maximum tensile stress. This is because 
increasing ao reduces (Cyyyy — Cyy^^) for angles near 0 = 45° (Figure 4a, black solid 
line). Superimposed on this is a compressive stress around 0 = 0° and 90° that emerges 
as a consequence of shear strain rate coupled to normal stress via the CyyxY viscosity 
(Figure 4a, gray solid line). The product CyyxY^^XY negative for all 0, motivating 
us to name this coupling “shear-strain-induced compression.” This non-trivial result 
comes from the fact that the stress-induced softening occurs in the tensile (ts) direction, 
as schematically illustrated in Takei and Katz [2013] Figure 7b. 

The effect of increasing 0o on Ty'^y is shown down columns in Figure 3. Similar to ao, 
00 couples the shear strain rate to compressive normal stress for angles near 0 = 0 ° and 
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Figure 3: Normal stress (tension positive) and shear stress T^y produced by the 
base-state simple shear flow as functions of angle 6 between the shear plane and the X-axis 
(Fig. 2a). We use = 5/3 here and throughout [Takei and Katz, 2013]. Each panel is 
computed with a different set (ao,/?o) as labelled above and to the right of the panels. 
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Figure 4: Viscosity components from equation ( 8 ) as functions of angle 6 between the 
shear plane and the X-axis. The curves are computed for the anisotropy parameters (a) 
oq = 1 and /3o = 0, and (b) oq = 0 and /3o = 1- The thin black lines represent each 
component for the isotropic case (uq = /3o = 0 ). 

90° via the CyY^y viscosity. In contrast to ao, however, Po strengthens the aggregate 
in the ri direction and increases the normal stress amplitude near 6 = 135° (Figures 3 
and 4b). As a result, the sign-change of normal stress caused by /do occurs in a limited 
range of 0 < 90° and 6 ^ > 0° (Fig. 3, bottom-left panel). 

The shear stress curves in Figure 3 (dashed lines) also change with increasing oq 
and/or (3o. For zero anisotropy in panel (a), the shear stress follows cos 20, as expected 
for coordinate rotation only. Anisotropy does not change the mean of T^y(0), as re¬ 
quired by symmetry of the stress tensor. Increasing (or /do) has an overall weakening 
(or strengthening) effect, changing only the amplitude of Tyy. This is in contrast to 
the effect of anisotropy on normal stress, which has a strong dependence on angle 6. 

4.2 Growth of porosity perturbations 

For simplicity in this linearised analysis, we consider the case of liquid viscosity rji = 0, 
giving a non-dimensional compaction length R ^ oo. In this limit, liquid segregation 
over any length scale occurs at vanishingly small pressure gradients. Therefore, the 
pressure gradient terms in equations ( 6 c) are negligible. Pressure then drops out of the 
problem and we no longer need to solve equation ( 6 b) (which has become indetermi¬ 
nate anyway!). This is equivalent to considering only the subset of perturbations with 
wavelengths much smaller than the dimensional compaction length [e.g. Katz et al., 
2006]. 

Using the initial wavenumber vector K = (0, K), the porosity perturbation at time 
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t is 

0i(X, t) = exp [iK ■ {X - V^^k) + St] , (9) 

which accounts for rotation of the wave-vector due to advection by the base-state flow 
[Spiegelman, 2003]. At t = 0, by choice of the coordinate system, perturbations are 
uniform in the X direction. Therefore, partial derivatives of the first-order quantities 
with respect to X are zero. Using = —e^y, = 0, and VPi = 0, equations (6) 
become 
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Equations (10b) and (10c) are obtained after an integration in the Y direction; boundary 
conditions are not needed because the domain is infinite and the first-order fields are 
periodic. 

The right-hand side of equations (10b) and (10c) represent Tyy and Txy, respec¬ 
tively. Since pressure gradients are negligible, these stresses must be spatially uniform 
for the system to be in balance. Therefore the first-order product of viscosity and 
strain rate must sum to zero; viscosity reduction associated with porosity perturba¬ 
tions (within the first square brackets in the right hand side) is compensated by the 
strain rate perturbations (within the second square brackets). 

To facilitate the physical interpretation of equations (10b) and (10c), these equations 
are re-expressed as 
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Equation (11) relates the strain-rate response of the system to the forcing stresses 
defined by (12). From equation (11), normal strain rate in the Y direction Cyy (the 
component that is most relevant to the perturbation growth) can be expressed as 
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The forcing stresses, Tyy and rj/y, are not externally applied (like those causing simple 
shear), nor are they the hrst-order stress perturbations Tyy and rj^y (these are both 
equal to zero). Instead, they are equivalent stresses that are created internally as a con¬ 
sequence of the base-state flow acting on the viscosity change associated with porosity 
perturbations. Moreover, under dynamic anisotropy, these forcing terms also depend on 
the strain-rate perturbations and hence equation (13) does not always give an explicit 
solution for Cyy. Nonetheless, equation (13) enables us to separate the mechanics into 
two simpler parts: the forcing, Tyy and r^y, and the compliance, Cyyyy and C'yyxv’ 
where the latter represents the system response to forcing with unit amplitude. This de¬ 
composition is helpful to understand the detailed (and rather complicated) mechanisms 
of the different models considered here. 


4.3 General solution 


Equations (10) are solved here to obtain an explicit expression for s for the full model of 
dynamic anisotropy. The hrst-order viscosity tensor is written in terms of the porosity 
and anisotropy perturbations 0i, ai, (3i, and ©i. The anisotropy perturbations are 
then expressed in terms of the porosity and strain-rate perturbations 0i, Cyy, and 
e^y. These calculations are sketched in Appendix B. The components of the hrst- 
order viscosity tensor are then substituted into equations (10b) and (10c), which are 
manipulated to solve for Cyy and e^y as functions of 0i. The normal strain rate 
Cyy obtained by this approach is substituted into (10a) to give an expression for the 
growth-rate of perturbations. 
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/^(O) J) I y 
^'YYyy^p '' 


YY h' XX j ' '^YYXY^q y XY P'XXj 


with the compliances dehned by (14) and dynamic factors 

1 “ P^xlcxy 

= ( 1 ~ Q^XXXy/^YYXY ] 


Vp = 


'\cYXY ] / 


where 


D = 1 + q 


-.( 0 ) 


^( 0 ) 


VO) 


VO) 


TW) TM ITVo) ;:Vo) 

^XYXY^YYVY ^YYXY^XYYY 


+ p 


VO) 


VO) 


^( 0 ) 


VO) 


TVo) y;(oi yTVo] 

^XYXY^YYYY ^YYXY^XYYY 


(15) 


(16a) 

(16b) 


(17) 


The constant coefficients p and q that express the sensitivity of the growth rate to 
dynamic anisotropy perturbations are 
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where a and ( are defined by equations (34) and (36), respectively. We do not attempt 
to physically interpret the detailed form of p and q. It is important to note, however, 
that for the static anisotropy model, p and q are zero, and Vp = Vq = 1. 

5 Physical interpretation in various limits 

The growth rate in eqn. (15) is a general result for the full model presented in §2 
above (with the sole assumption of i? —)■ oo). To build up a physical understanding 
of this equation, we return to the simpler case of static anisotropy, which includes 
the simplest case of the Newtonian, isotropic model. The static anisotropy model has 
previously been studied by Takei and Katz [2013] using linearised analysis. However, the 
mathematical complexity of their results precluded a detailed mechanical interpretation. 
A reconsideration using the perturbation-oriented coordinate system enables a physical 
understanding of the instability mechanism and the rheological control on the dominant 
band angle. These are needed to understand the more complicated, dynamic model. 
To facilitate this (in §5.1-5.2), we make the simplifying assumption that (3 = 0 — that 
there is no contiguity increase in the ri-direction. The effect of non-zero (3 is discussed 
in section 5.3, where we show that its role is minor compared to that of a, the contiguity 
decrease in the rs-direction. 

5.1 Static anisotropy 

When «! = /9i = ©1 = 0, the mechanical equilibrium conditions (10b) and (10c) are 
written as 


Cyyyy^yy + (19a) 

/^(O) a(1 ) pd) _ _(0) ('19b') 

XYYY^YY ' ^^XYXY^XY ~ 

/ 

Comparison with equation (11) shows that the forcing stresses are given by Tyy = 
A^iTyy and t^y = ^4>iTxy- These forcing stresses are caused by the base-state tensile 
and shear stresses acting on the porosity perturbation by way of porosity weakening 
rheology (A > 0) as depicted in Figs. 2b, c. In this simple model, Ty^ and r^y are 
given in terms of the porosity perturbation 0i, and hence eqn. (13) provides an explicit 
solution for Cyy. 

It is evident from (10a) that the normal strain rate Cyy causes an increase in the 
amplitude of porosity perturbations; the shear strain rate e^y does not cause the poros¬ 
ity to change. When the viscosity is anisotropic, both the normal and the shear stress 
drive Cyy and hence contribute to perturbation growth. This does not occur under 
isotropic viscosity. 
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5.1.1 Instability mechanism in the isotropic system 

For an isotropic aggregate (oq = /^o = 0), Cyyxy — ^xyyy = 0 in equations (11) and 
(19), and the compliance CyyxY couples shear stress to normal strain rate is zero. 
In this case, eyy is driven only by the base-state normal stress Tyy. The growth rate 
s is given by 

^( 0 ) 

s = (1 — 0 o)A—— isotropic model. ( 20 ) 

^YYYY 

Perturbations are unstable under tensile stress {Tyy > 0) normal to the perturbation 
wavefronts, stable under compressive normal stress (r^y < 0 ), and unaffected by shear 
stress T^y. We therefore term this the tensile-stress-induced instability, or tensile insta¬ 
bility. When band angle 6 relative to the simple shear flow is 45°, the tensile stress Tyy 
attains its maximum (Figure 3a) and hence the growth rate s is also at a maximum, as 
shown in Figure 5a and by Spiegelman [2003]. The occurrence of the tensile instability 
in a porosity-weakening, two-phase aggregate was hrst predicted by Stevenson [1989]. 

5.1.2 Two instability mechanisms in the anisotropic system 

For an anisotropic aggregate (oq > 0 and/or > 0), there is a coupling between shear 
and normal components via Cyyy.y = C^yyy ^ 0. In this case, Cyy is forced by both 
normal stress across perturbations and shear stress along perturbations. The growth 
rate is 

s = (1 — (j)o)X (CyyyyTyy + C'yyxy'^xy) static anisotropy model, (21) 

using the compliances given by equations (14). The hrst term on the right-hand side of 
(21) represents the tensile instability, generalised to the anisotropic aggregate. The 2nd 
term represents a shear-stress-induced instability that does not occur in the isotropic 
system. The total growth rate s versus band angle 6 is plotted in the top row of panels 
of Figure 5 for various anisotropy amplitudes (thick lines). Consistent with previous 
work, as is increased, the single growth rate peak splits into two peaks at low and 
high angles to the shear plane (Fig. 5c). Because the lower-angle peak dominates the 
higher angle peak after a hnite time [Katz et al., 2006, Takei and Katz, 2013], this result 
means a signihcant lowering of the dominant band angle by the viscous anisotropy — 
if the magnitude of anisotropy a is sufficiently close to saturation (a ~ 2 ). 

5.1.3 How viscous anisotropy causes lowering of band angle 

Although the effect of viscous anisotropy is evident from the total growth rate shown 
in the top panels of Figure 5, it is not immediately obvious why the dominant band 
angle is lowered by viscous anisotropy. The physical mechanism can be understood by 
considering the tensile and shear components of the growth rate independently ( 1 st and 
2nd terms of (21), respectively). In the top row of Figure 5, these two growth rates 
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ao = 0 


ao — 1 


ao — 2 




0, degrees 0, degrees 0, degrees 


Figure 5: Characteristics of the static anisotropy model as a function of the angle between 
porosity perturbations and the shear plane. Each column is for a different value of ao- In 
all cases, /3o = 0, A = 27, and = 0.05. Top row. Growth rate s of perturbations c/)! 
from eqn. (21). The heavy line represents the total growth rate; light lines show the growth 
rate decomposed into two parts: the contribution from the tensile instability (1st term 
of eqn. (21)) and that from the shear instability (2nd term). Bottom row. Compliances 

Cyyyy ^YYXY ™ ( 21 ). 
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are plotted separately for various values of ao (thin solid curve for tensile instability; 
thin dashed curve for shear instability). Comparison of panels (a), (b), and (c) reveals 
that the peak split occurs through (z) stabilisation of the tensile instability and (n) 
emergence of the shear instability with increasing magnitude of anisotropy a. We 
consider each of these in turn. 

To understand why viscous anisotropy stabilises the tensile instability, we return to 
the systematics of the base-state stress (section 4.1). Comparison of the three columns 
of Fig. 3 shows that as ao increases, the tensile stress Tyy decreases in amplitude and 
becomes compressive at all angles. With Tyy < 0, the hrst term in equation (21) is 
always less than or equal to zero, and hence stable. 

To understand why viscous anisotropy destabilises the shear mechanism, we consider 
the coupling between the shear stress that drives the instability and the normal strain- 
rate that is responsible for its growth. As shown by equation (13) with Tyy = A0iryy 
and T^y = X(j)iT^Y, the shear stress r^y is coupled to the normal strain rate eyy 
via Cyyxy- The angular dependence of is shown by dashed curves in the 

bottom-row panels of Figure 5. If CyyxyTxY is positive, then s is positive (or Cyy is 
in phase with 0i) and the shear mechanism contributes to unstable growth of porosity 
perturbations. In fact, this product is positive (or zero) for all 6, enabling us to name 
this coupling “shear stress-induced expansion.” This non-trivial result comes from 
the assumed microstructural behaviour: that stress-induced softening occurs in the 
tensile (ts) direction, as illustrated in Takei and Katz [2013] Figure 7a. So the porosity 
perturbation grows because of the shear mechanism, for which the low angle is favorable. 

5.1.4 Summary of static anisotropy model 

As a recap and summary, note that under isotropic viscosity, the growth of bands at 45° 
to the shear plane is caused by a tensile instability [Stevenson, 1989, Spiegelman, 2003]. 
In contrast, under anisotropic viscosity, the peak growth rate of bands is controlled by a 
distinct shear instability. Although the peak growth rate of the shear instability occurs 
at 6 < 15°, stabilisation at these low angles by the tensile mechanism acts to give a 
maximum in the combined growth rate at 6^ = 15°. 

The comparison between isotropic and anisotropic systems developed above is sum¬ 
marised in the hrst two rows of Table 1. The tensile instability is separated into 
porosity-weakening A, which is fundamental to all models, and the tensile stress across 
bands Tyy, which affects both isotropic and anisotropic cases. A shift of Tyy to more 
negative, compressive values (represented by A) stabilises the tensile instability. In 
contrast, the difference in shear instability can be simply shown by its existence or non¬ 
existence (O or ^)- It is the leading-order terms Tyy and C^yxy that are responsible 
for these differences. 

Katz et al. [2006] extended the analysis of isotropic viscosity to include a power- 
law dependence of viscosity on strain rate (or, equivalently, on stress). They showed 
that strain-rate weakening viscosity leads to lowering of band angle. In the discussion 
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Table 1: Summary of band formation models. Q- exists; does not exist; A: is modified. 


Model 


Isotropic Newtonian 
Static anisotropy 
Dynamic anisotropy 
Isotropic power-law 


Cyyyy'^yy (tensile) 

A > 0 Tyy > 0 additional 
factor 

o 

O A 

O A o 

o o o 


^YYXY "^XY (shear) 


o 

o 


Dominant angle 


45° 
~ 23° 
~ 10 ° 
~ 20 ° 


section, we compare the angle-lowering mechanism of viscous anisotropy to that of the 
power-law viscosity. This is enabled by a reanalysis of the power-law model using the 
rotated coordinate system. 


5.2 Dynamic anisotropy 

We return to the full expression for the growth rate of bands, eqn. (15), to develop 
a physical understanding of why dynamic anisotropy lowers band angles, as observed 
in the numerical solutions (Fig. 1). To do so we take ao = 1 and again make the 
simplifying assumption that (3 = 0 (though see §5.3, below). 

The perturbations in 0 and a under dynamic anisotropy are obtained by linearisa¬ 
tion of equations (3) and (4) with respect to the stress perturbation . The expansion 
is conducted around the base-state values 0o,ao- In Appendix B, we show that the 
sensitivity of a to variations in deviatoric stress is given by the parameter 


da 


a = 


dAr 


= -sech^ 

Ar(®) "^sat 


2(Ar(°^ - Toffset) 
"^sat 


( 22 ) 


This parameter allows us to write cki = d, and hence to see that static anisotropy 
corresponds to the case where d = 0. The situation for 0 is more complicated because 
there is no single parameter that controls its sensitivity to deviatoric stress; varia¬ 
tions of 0 can either be fully considered or fully neglected. Fortunately, numerical and 
analytical results show that these variations (0i) play an insignificant role in the under¬ 
standing of band angles, and hence we need consider only the magnitude of anisotropy 
a. This is achieved by looking at the dependence of key quantities (especially s) on d. 

The growth-rate of porosity perturbations s is shown in panel (a) of Figure 6 for 
ckq = 1 and for values of d ranging from zero to two. Although d = 0 does not 
exclude linearised variations in 0, comparison with the s-curve in Fig. 5b conhrms 
that variations in 0 are ineffectual; with d = 0 the expected band angle is 45°. For 
increasing d, the growth-rate peak again splits into peaks at low and high angles. It is 
important to note that the mean value of a (cto) is not changed in this exercise — only 
the amplitude of variations about that mean. Consistent with the numerical results of 
Figure 1, dynamic variations in the magnitude of anisotropy can sharply reduce band 
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Figure 6: Characteristics of the dynamic anisotropy model for various values of d as a 
function of the angle 9 between porosity perturbations and the shear plane. In all panels, 
ao = 1) /3o = 0, A = 27, and (/>o = 0.05. (a) Full growth-rate s from eqn. (15). (b) 
Growth-rate from eqn. (15) split into the tensile-instability term (solid lines) and the shear- 
instability term (dashed lines), (c) Band-normal forcing stress (solid lines; eqn. (12a)) 
and band-parallel perturbation stress (dashed lines; eqn. (32)). Both are divided by 
Calculation details are in Appendix B. (d) and (e) The anisotropy perturbations ai and 
01 divided by Aijii, calculated according to eqn. (35). 
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angles, even at moderate ao for which the static anisotropy model predicts a high band 
angle (45°). 

Figure 6b breaks the full growth rate into two parts, each associated with one of 
the terms of equation (15). Dashed lines, representing the shear instability, are almost 
unaffected by a. In contrast, the tensile instability is strongly stabilised with increasing 
a. This stabilisation causes the peak of the full growth rate in panel (a) to split into 
low- and high-angle peaks. To understand why dynamic anisotropy promotes low band 
angles, it is therefore sufficient to understand why it stabilises the tensile instability. 

The tensile instability is driven by Tyy, as discussed in §5.1. This represents the 
normal stress (tension positive) that arises when viscosity perturbations interact with 
the base-state strain rate. The detailed form of the forcing stress for the dynamic 
anisotropy model is given in equation (31). Figure 6c shows that the forcing normal 
stress Tyy varies signihcantly with a (whereas the forcing shear stress, not shown, is 
almost unaffected by a). The system compliances, which are leading-order quantities, 
are not affected by dynamic anisotropy. Therefore, it is the variation of r-py that is 
responsible for stabilisation of the tensile instability under dynamic anisotropy. 

To develop a physical understanding of the detailed dependence of Tyy on 01; cq 
and ©1 (eqn. (31a)), focus attention on 6^ = 45°, as this is the dominant band angle 
when a = 0. For bands at 45°, solid curves in panel (c) of Figure 6 show that the forcing 
stress goes from a positive perturbation (in phase with 0i) to a negative perturbation 
(anti-phased with 0i) with increasing a — hence the forcing stress in the high- 
porosity bands goes from tensile to compressive. This change is due to an increase in the 
magnitude of anisotropy perturbation ai = aAT^^\ shown in panel (d). Since t^y = 
Tyy = 0, the deviatoric stress perturbation is entirely due to the band-parallel 

normal stress perturbation (according to eqn. (32)), which is shown by dashed 
curves in panel (c). Because < 0, > 0 signihes a magnitude reduction of txx 

in the high-porosity bands; the largest change occurs for bands at 6^ = 45°. As sensitivity 
a to deviatoric stress increases, ai becomes more negative (panel (d)). Negative values 
of (anti-phased with 0i) mean high-porosity bands have lower deviatoric stress and 
weaker anisotropy than the low-porosity, inter-band regions. This is consistent with 
numerical results in Fig. Id. 

Figure 5d-f shows that ao increases the normal compliance Cyyyy at angles between 
zero and 90°. A negative perturbation to ao therefore makes the high-porosity bands 
in this range of angles less compliant to tensile stress and the low-porosity inter-bands 
more compliant. Overall, then, the perturbation in anisotropy amplitude ai tends to 
cancel the direct effect of the porosity perturbation 0i on the normal compliance, and 
hence ai works to stabilise the tensile instability. 

The comparison between the static and dynamic anisotropy models developed in this 
section is summarised in Table 1. These two models are identical at leading order but 
different at hrst order. Therefore, stabilisation of the tensile mechanism due to more 
compressive base-state stress (ryy) and destabilisation of the shear mechanism due 
to shear stress-induced expansion {Cyyxy) occur in both the static and the dynamic 
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anisotropy model. These two cases differ, however, in that further stabilisation of the 
tensile mechanism occurs due to the dynamic variation of anisotropy magnitude (ai). 
This effect hardens the band regions and weakens the inter-band regions under dynamic 
anisotropy. This additional factor (Q in Table 1) significantly lowers the band angle. 

It is interesting to note that dynamic perturbations to the angle of anisotropy ©i are 
not an important control on band angle. Figure 6e shows that they are not affected by 
a. More importantly, ©i is always zero for bands orientated at 6^ = 45°. This indicates 
that the stabilisation of the tensile instability and the lowering of band angle under 
dynamic anisotropy cannot be attributed to ©i. In numerical simulations (Fig. le), 
the variations of © do not contribute to the lowering of band angle that is observed in 
Fig Ic, though they are well-explained by the stability analysis at 6^ ~ 10° (red dashed 
line). 


5.3 The effect of contiguity increase in the ri-direction 

Until now, we have neglected {3 and focused on the effects of a, which quantifies conti¬ 
guity decrease in the direction of maximum tension. Non-zero a represents a weakening 
in the Ts direction that {i) reduces the magnitude of tensile stress and leads to (m) shear 
strain-induced compression and (in) shear stress-induced extension [Takei and Katz, 
2013]. We have shown that the tensile mechanism is stabilised around 6 = 45° by the 
first of these and is stabilised around 6 = 0° and 90° by the second; we have also shown 
that the shear mechanism is destabilised by {Hi). The parameter jS quantifies the con¬ 
tiguity increase in the direction of maximum compression. Even if a is zero, a non-zero 
(3 creates viscous anisotropy (see eqn. (2b)), causing the couplings {ii) and {in). How¬ 
ever, Fig. 4 shows that (3o does not cause the weakening (z). It is this weakening, by a 
only, that is responsible for splitting the growth-rate peak in both static and dynamic 
models. On this basis, we expect the effect of (3o to be small. This is indeed the case: 
as shown below, f3 alone does not cause a lowering of band angle, but it can affect the 
lowering by a. 

The left column of Figure 7 illustrates the effect of [3 under static anisotropy. 
Panel (a) shows that under static anisotropy, s is split into high and low angle peaks for 
any value of /Sq = 0 when Oq = 2 (solid curves), whereas it is peaked at 45° for any value 
of (3q when ao = 0 {Pq = 2 shown by dash-dot curve). For ao = 2, increasing (3o causes 
a modest shift to more compressive Tyy at 0 ~ 0° and ~ 90° and a modest increase in 
the amplitude of shear stress r^y (Fig. 3, right column). Therefore, as Figure 7b shows, 
Pq causes stabilisation of the tensile instability and destabilisation of shear instability 
in equal measure. These two effects compensate each other and the solid growth-rate 
curves in Figure 7a are thus all very similar to that for f3o = 0. 

The right column of Figure 7 shows how 13 affects dynamic anisotropy. Panel (c) 
shows that for s in the dynamic anisotropy model, a two-peaked growth rate occurs 
for ao = d = 1 and (3q = $ = 0 (light gray curve) but does not for ao = d = 0 and 
/do = /3 = 1 (dash-dot curve). In the former case of non-zero a with a double peak, 
increasing (3q enhances the stabilisation of the tensile mechanism at 45° and deepens the 
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Pq = 2, oq = 0 
/3o = 0, oq = 2 
' 1^0 = li oo = 2 
• Pq = 2, aa = 2 


Static 


Po 

Po 


P = 1, oo = a = 0 
P = 0, oo = a = 1 

/3 = l/2,ao = a = l 


Dynamic 


Po = p = l, Oo = a = 1 



Figure 7: The effect of /3 on the growth rate of porosity perturbations for static (left 
column) and dynamic (right column) anisotropy, (a) Full growth rate s from eqn. (21) for 
ao = 2 and various values of /3o- A single curve for oq = 0 and /3o = 2 is also shown. In 
panel (b), the cases with cto = 2 are decomposed into tensile and shear parts. Line greyscale 
has the same meaning as in panel (a); there is no curve for oq = 0. (c) Full growth rate s 
from eqn. (15) for oq = « = 1 and various values of /3o = /3. A single curve for oq = « = 0 
and /3o = /3 = 1 is also plotted. In panel (d), cases with = a = 1 are decomposed into 
tensile and shear parts, (e) The band-parallel, normal stress perturbation 
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valley between low-angle and high-angle peaks of s (Fig. 7c-d). This occurs because Txx 
is enhanced by the overall strengthening effect of Po (Fig. 7e). The very low band angles 
that emerge in the numerical simulation with dynamic anisotropy and = 1 (Fig. 1) 
are therefore a consequence of both the dynamic effect of a and the enhancement by (3. 

6 Summary and discussion 

We have developed and analysed a model of coupled magma/mantle dynamics with 
anisotropic viscosity. The anisotropy is controlled by the orientation of principal stresses 
and the amount of deviatoric stress. The model presented here introduces small mod- 
ihcations on that of Takei and Katz [2013]; in particular, the parameter f3 models an 
increase in contiguity of grains in the direction of maximum compressive stress and the 
parameter Toffset allows for a hner control on the magnitude of anisotropy and its sen¬ 
sitivity to stress (for = const.). This description of viscous anisotropy is physically 
consistent with experiments and relatively simple, so its analysis should clarify the me¬ 
chanics of rocks for which the assumptions hold. Existing experimental data, however, 
are not enough to quantitatively constrain all parameter values. The parameter studies 
performed here aim to understand the underlying physics. 

It is known from previous theoretical work that anisotropic viscosity lowers the 
angle of emergent, high-porosity bands. Numerical solutions (Fig. 1) compare uniform 
anisotropy imposed a priori with anisotropy that varies according to local conditions of 
stress. They show that dynamic anisotropy leads to lowering of band angle as compared 
with uniform anisotropy, where the mean magnitude and angle from the dynamic case 
are used in the static case. Moreover, dynamic anisotropy produces low-angle bands 
even when its mean values wouldn’t do so if applied uniformly and held constant. The 
physical reasons for this have not previously been clear. Indeed, the question of why 
anisotropic viscosity lowers band angle at all has not previously been addressed. 

Static viscous anisotropy, in which viscous resistance to extension in the most tensile 
direction is decreased, predicts low-angles of high-porosity bands for two reasons: (a) it 
suppresses the mode of instability in which tension causes extension across high-porosity 
bands; (6) it creates a mode of instability in which shear stress causes extension across 
high-porosity bands. The tensile instability has a peak perturbation growth rate in 
the maximum tensile direction {9 = 45°). When this instability is suppressed by static 
anisotropy, the peak growth rate shifts to the smaller angles that are favoured by the 
emergent shear instability. And although the growth of the lowest angle bands are 
enhanced by the shear instability, perturbations parallel to the shear plane {9 = 0°) 
are stable because of the compressive stress created by the base-state flow. Therefore, 
a low but hnite angle of high-porosity bands is predicted by this model. Allowing for 
an increase in contiguity and viscosity in the direction of maximum compression has 
counter-balancing effects that leave predicted band angles almost unchanged. 

Dynamic viscous anisotropy, in which the anisotropy parameters are allowed to vary 
with the local orientation and magnitude of deviatoric stress, tends to further lower 
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band angles. It does so because it suppresses the tensile instability around 6 = 45° 
via the following dynamic effect. Lower deviatoric stress in viscously weak bands gives 
lower anisotropy there, which makes them less compliant to tensile stress across them. 
Enhanced anisotropy in the interleaved, lower-porosity regions makes those regions 
more compliant. This effect over-compensates the compliance variations directly due 
to porosity weakening; it favours melt segregation from the bands into the inter-bands. 
Allowing for an increase in contiguity and viscosity in the direction of maximum com- 
pressional stress increases the contrast in band-parallel compressional stress (and devi¬ 
atoric stress) between bands and inter-bands. This enhances the contrast in anisotropy 
and further suppresses the tensile instability. Dynamic anisotropy makes almost no 
modihcation to the shear instability. 

The additional effects of dynamic anisotropy and the anisotropic increase of conti¬ 
guity are important because they make more robust the prediction of low band angles. 
Under static anisotropy, the mean magnitude of anisotropy must be quite high to pro¬ 
duce low-angle bands; moderate levels are insufficient. In contrast, under dynamic 
anisotropy with contiguity-increase in the direction of maximum compression, moder¬ 
ate levels of mean anisotropy efficiently produce low-angle bands. This helps to support 
the hypothesis that low-angle bands in experiments are due to anisotropic viscosity be¬ 
cause it expands the parameter space in which the theoretical predictions should hold. 

These conclusions were reached by use of stability analysis in a coordinate system 
that is rotated with respect to the plane of simple shear; in particular, the coordinate 
system is aligned with the wavefronts of the harmonic perturbations. This rotation 
leads to simpler expressions for the growth rate of perturbations: the tensile and shear 
modes appear as distinct terms that are amenable to physical understanding. For this 
reason, our analysis represents a framework in which to test and understand the family 
of rheologies that potentially produce low-angle bands in shearing flows. This includes 
variants of isotropic and anisotropic viscosity, but also potentially of dilatational gran¬ 
ular rheology, damage, or composite rheologies [e.g. Rudge and Bercovici, 2015]. 

An application of the rotated coordinate system to the isotropic power-law creep 
model with stress exponent n [Katz et al., 2006] is presented in Appendix A. As with 
all other models considered here, porosity perturbations reduce viscosity in the bands, 
resulting in the enhancement of the normal and shear strain rates, eyy ^md e^y. In this 
model, however, the enhanced strain rates further reduce the viscosity, which feeds back 
to further enhance the strain rates. The importance of this non-Newtonian feedback 
relative to the porosity-weakening feedback is roughly approximated by n — 1. Both 
normal and shear strain rates contribute to the non-Newtonian feedback; the relative 
importance of the shear component increases with increasing r^. Therefore, if n and 
are sufficiently large, shear strain rate is the key weakening factor and the growth rate 
reaches a maximum at a substantially lowered angle. However, in contrast to anisotropic 
viscosity, strain-rate weakening viscosity does not give rise to the shear instability — it 
merely lowers the most favorable angle for the tensile instability (comparison in Tab. 1). 
Although the details differ, both models predict an important role for shear stress in 
the lowering mechanism; both predict a low but hnite angle with localised shear strain 


24 


in the higher-porosity bands. 

The model of viscous anisotropy used here seems promising as an explanation for 
laboratory experiments on deformation of partially molten rocks. Although its detailed 
form must be considered tentative, we are not aware of another theory that reproduces 
the low-angle bands found in experiments [Holtzman et al., 2003] while respecting the 
measured stress-dependence of creep viscosity [King et al., 2010]. Furthermore, radially 
inward migration of magma in experiments employing torsional deformation [King et al., 
2011, Qi et ah, 2015] may be direct evidence of base-state segregation, a feature that 
arises naturally from viscous anisotropy [Takei and Katz, 2013] but may be impossible 
to reconcile with isotropic viscosity. Although the present study focuses on the angle of 
bands, the growth rate of bands is also affected by static and dynamic anisotropy; the 
growth rate is lowered by static anisotropy and further lowered by dynamic anisotropy 
(Figs. 5 and 6). This can be also discerned in the different total strain and different 
ranges of porosity in Figs, la and lb. Therefore, a quantitative comparison between 
the measured and predicted growth rate becomes important for further rehning and 
testing the theory. 

In the present theory, a, /9, and 0 are assumed to depend on stress, based on 
the experimental results by Daines and Kohlstedt [1997] and Takei [2010]. Although 
this assumption is considered to be valid at small strain, possible evolution of these 
parameters with increasing strain has to be investigated to model the system at large 
strains. Indeed, for more than 200% strain under simple shear, Zimmerman et al. [1999] 
observed that the long axis of melt pockets is predominantly oriented at an angle of 
20° from ri; this is difficult to explain by stress alone. It should be noted, however, 
that microstructural analysis in laboratory studies has been performed in terms of 
shape and orientation of melt pockets; an analysis in terms of observed contiguity is 
more appropriate for comparison with and incorporation into the model. Numerical 
simulation using dynamic anisotropy and an empirically justihed evolution equation for 
contiguity will be important in future work. 
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A Power-law creep model by Katz et al. (2006) 


The model of band-formation under power-law viscosity by Katz et al. [2006] is formu¬ 
lated by equations (1) and the viscous constitutive relations 


iji kl^XX YY XY 
XX/ r, + | r,-| 0 \ 

Cijki = r/(0, eii) X YY \ ■ | 0 h 

xrV ■ • 1 / 


(23) 


where only 6 of the 16 components of the two-dimensional version are shown due to the 
symmetry of Cijki- The normalised shear viscosity ri{(j), en) depends on porosity and 
the second invariant of the strain-rate tensor, en = as 

1 —n 

r]{(j), en) = exp [-A((/) - (j)o)/n] ejf (24) 


[Katz et ah, 2006, Takei and Holtzman, 2009b]. Equation (24) represents a power- 
law viscosity that, to represent deformation by dislocation creep, has an exponent 
n ^ 3.5 [e.g. Karato and Wu, 1993] (n = 1 corresponds to Newtonian viscosity). Under 
dislocation creep, the strain rate is highly sensitive to the stress because dislocation 
velocity and density both increase with increasing stress. Hence the model of Katz 
et al. [2006] incorporates strain-rate weakening in addition to the porosity weakening. 

Katz et al. [2006] demonstrated that strain-rate weakening viscosity works to lower 
the band angle. Although the mechanism of this lowering is briefly discussed in their 
paper, further analysis of their model using the perturbation-oriented coordinate system 
is helpful to understand their explanation and to compare it with the mechanism of 
viscous anisotropy. For consistency with the foregoing development, = 0 is assumed 
here. We can expand (24) into base-state and perturbation terms as 
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1 —n 

where r]o = we have used e® = 1/2 and e^^x ~ 0- Combining 

(23) and (25) with stress balance (10b) and (10c), we obtain 
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(26a) 

(26b) 


This formulation is not the most amenable to inversion for the strain-rate perturbations, 
but it allows for a clear comparison with equations (11) and (12). We have moved terms 
to the right-hand side that can be considered to comprise the forcing stresses rpy and 
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9, degrees 


Figure 8: Characteristics of the power-law isotropic viscosity model as a function of the 
angle between porosity perturbations and the shear plane. In both panels, n = 5. (a) 
Normal strain rate perturbation, which shows a double peak. (b) Forcing normal stress 
Tyy due to porosity weakening (solid curve), strain-rate weakening associated with eyy 

(dashed curve), and strain-rate weakening associated with e^y (dotted curve). The three 
curves in panel (b) sum to the curve in panel (a). 


Txy- Two points are evident: First, the forcing stresses retain the term representing 
base-state stress operating on porosity perturbations. Second, there are new terms that 
cross-couple the equations (26). 

The cross-coupling terms in (26) arise because normal eyy and shear cxy compo¬ 
nents both affect e// and hence modify the viscosity (by way of an increase in dislocation 
density). Two feedback mechanisms are thus at work, causing growth of porosity per¬ 
turbations. The hrst of these is a direct effect: when A > 0, high-porosity bands are 
weaker by virtue of their higher porosity. The second is indirect: porosity-weakened 
bands have a larger strain-rate that, when n > 1, further weakens them through the 
non-linear viscosity. The relative importance of the 2nd mechanism to the 1st one 
increases with increasing n — 1. 

Solving equations (26) for Cyy and e^y and using equations (10a) and (23), the 
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growth-rate is 
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Takei and Holtzman [2009b] obtained the identical result for s and showed that a single 
peak splits into two at large n. Here, to understand the mechanism of the split, the 
forcing stress associated with tension, r-py, is plotted in Figure 8b for each of the three 
terms in the RHS of (26a). Although the tensile forcing stress due to the porosity and 
normal-strain rate perturbations are maximum at 6* = 45° (solid and dashed curves), 
that due to the enhanced shear strain rate e^y has peaks at 0 ~ 10° and 80° (dotted 
curve). The sum of these three curves determines the prohle of Cyy in panel (a) and 
hence determines the growth-rate. Panel (b) conhrms that the weakening of viscosity 
by the enhanced shear-strain rate e^y is the main cause of the peak split of the growth 
rate. 


B Calculation of the dynamic-anisotropy growth rate 


To solve the hrst-order equations of force balance (10), we need an expansion of the 
viscosity tensor into its base-state and perturbation components. Equation (2b) gives 
Cijki as a function of 0, a, 13, and 0. Under the dynamic anisotropy model, it is 
necessary to account for non-zero perturbations CKy /9i, and ©i. In that case, and 
C-jh are calculated as 


— Cijki{4>o, ao, (3o, ©o) 
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(28b) 


From the equation for anisotropy magnitude (4a), 


cto = 1 + tanh 


/ 2ArW-2roffset \ 

\ "^sat / 


(29) 


where = 2(1 — (ao ~ /So)/4). Then from equation (28b), the stress perturbation 

is written as 
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Using (30), (2b), (7), and ©o = vr/4 + 6*, the mechanical eqnilibrinm conditions = 
Tyy = 0 from eqnations (10) are written as 

CyyyyCyy + 2CyYXY^^XY = ^ (“l 26* + 200©! COS 26) + ^ 

(31a) 

Cj^yyyCyy + 2uj|’yjf yc^y = A^iT^^y H-(«! COS 26* — 2ao©i sm 26*), (31b) 

where the right-hand sides of these eqnations are the dynamic-anisotropy version of 

f/’) { 

the forcing stresses r-yy and rj/y, respectively. The normal-stress perturbation in the 
X-direction is 


,(1) _^(0) •(!) , 

'XX — '-'XXYY^YY ^ ‘^'-'XXXY^XY '^Yl'xX 

— {ai sin 26* -|- 2 qo©i cos 26*). (32) 

Microstructural anisotropy is determined by deviatoric stress. From the total dif¬ 
ferentials of equations (3) and (4a), and from r^y = ^yy = 0, ai and ©i are related to 
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Then we use the expression (32) for Txx and equations (33) to obtain 
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where 


r^ = 


l + r0 


r0 


sin 26 


2Ar(°) 


a[T 


(0) 


r, 


(0)n 

yy) 


1 - 


(cTo ~ /^o) 


°vxx ^yyJ 


( 0 ) 


cos 26* sin 46* 


(36) 

Equations (35) can be substituted into the stress-balance equations (31) giving a system 
in which the only hrst-order quantities are 0i, e^y, and Cyy. 
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